clear all
numIter = 5000;
seqLen = 128;
period_labels = [7:15];
numTargets = 7;
ITI = 3200;

maxDetectPower = 0;
best_baseSeq = [];
best_guidefunction = [];
best_fullSeq = [];
best_detectPower = [];


for i=1:numIter
    fprintf('Iteration %d: ',i);
    % Generate sequence using a sinusoidal guide function with period=period_labels
    [baseSeq, guideFunction] = GenerateLongTermSequence(seqLen,period_labels);
    
    % Add zeros and targets
    fullLen = seqLen + numTargets;
    
    % Choose a subset of indx to be targets, avoiding the first and last 10
    targetIndices = randperm(fullLen-20)+10;
    targetIndices = sort(targetIndices(1:numTargets));
    
    % Make sure targets are spaced at least 3 trials apart and last target
    % is within the last 20 trials
    while (min(abs(targetIndices - [targetIndices(2:end) targetIndices(1)])) <= 3) || (abs(max(targetIndices) - fullLen) > 20)
        targetIndices = randperm(fullLen-20)+10;
        targetIndices = sort(targetIndices(1:numTargets));
    end
    
    % Add targets to baseSeq
    fullSeq = baseSeq-1;
    for thisTarget = targetIndices
        fullSeq = [fullSeq(1:thisTarget-1) -1 fullSeq(thisTarget:end)];
    end
    
    
    % Calculate detection power
    detectPower = evalSeq(fullSeq, ITI, 'hrf.txt', 10);
    
    if detectPower > maxDetectPower
        
        best_baseSeq = baseSeq;
        best_guidefunction = guideFunction;
        best_fullSeq = fullSeq;
        best_detectPower = detectPower;
        
        if detectPower > maxDetectPower
            maxDetectPower = detectPower;
        end
    end
    
    fprintf('maxDP = %f \n',maxDetectPower);
end


save longTermSeq